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The dynamics of particle transport under the influence of localised high energy anomalies (explo¬ 
sions) is a complicated phenomena dependent on many physical parameters of both the particle and 
the medium it resides in. Here we present a conceptual model that establishes simple scaling laws for 
particle dispersion in relation to the energy released in a blast, properties of the medium, physical 
properties of particles and their initial position away from a blast epicenter. These dependencies 
are validated against numerical simulations and we discuss predictions of the model which can be 
validated experimentally. Other applications and extensions to the framework are also considered. 
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I. INTRODUCTION 

Understanding transport properties of particle systems 
driven by strong energy fluxes is of significant importance 
to a number of fields of science and technology. Ex¬ 
amples include non-equilibrium statistical physics Una, 
astrophysical and geophysical phenomena m, multi¬ 
phase turbulent flows [B], inhomogeneous catalysis [7J, 
combustion [Sj and many others mm- An important 
requirement for these studies is the development of a 
rigorous framework which estimates system parameters 
( e.g. energy fluxes) from remote (or retrospective) obser¬ 
vations of particles following natural or anthropogenetic 
phenomena (a volcano eruption, meteorite impact, su¬ 
pernova event, blast, etc). Two revealing examples of 
this approach are the well-known pioneering studies of 
L.F. Richardson (estimation of parameters of turbulent 
flows from Lagrangian measurements mi) and E. Fermi 
(a remote estimation of a nuclear bomb yield from ‘tracer’ 
particle observations ED- 

There is a vast amount of literature devoted to the 
subject of particle transport (see P [BHS Eg and ref¬ 
erences therein). Modern computational models (often 
called models of Lagrangian transport) achieve an un¬ 
precedented level of fidelity by matching numerical pre¬ 
dictions with experimental observations 0H2Q2]. Un¬ 
fortunately, whilst these computational models are an 
important predictive tool for practical applications and 
validation studies, they are unable to provide analyti¬ 
cal insights into the fundamental transport mechanisms 
of these systems - simply because analytical predictions 
cannot be deduced numerically. Significant analytical 
progress in the understanding of transport phenomena 
in particle systems has been achieved by employing scal¬ 
ing and self-similarity frameworks mmm- This allows 
us to describe the dispersion process by means of power- 
law functions (scaling laws) relating to particle displace¬ 
ment and other parameters of the system - the exponent 
of these power-laws being predicated analytically. Using 
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this knowledge in conjunction with Lagrangian measure¬ 
ments, one can infer values of important system param¬ 
eters that would be challenging to recover by any other 
means. 

The presented results are in line with this approach. 
More specifically: we establish a scaling law for particle 
dispersion caused by a rapid (and localized) energy re¬ 
lease (explosion), and express it as a simple power-law 
relation between the particle displacement and the phys¬ 
ical parameters of the system (energy of the explosion, 
particle properties and their initial position, etc). We 
demonstrate that under a broad range of conditions the 
exponents of this scaling law can be deduced analytically 
by applying the ideas of self-similarity. We support our 
analytical predictions with numerical simulations. 

II. MODEL 

The following is a simplified conceptual model that al¬ 
lows us to derive a scaling law for particle displacement. 
We consider an infinite domain initially populated with 
particles whose density is much larger than the media 
density (we assume this media to be a gas with known 
properties). We restrict ourselves to the case when the 
density fraction of particles is relatively small, so particle- 
particle interactions can be neglected and thus consider 
the dispersion of a single tracer particle (the opposite 
limit of a localized energy release in a ‘crowded’ system 
of particles was analysed in Ref. [T9] Similar to other 
studies }20j . we assume that the dynamics of an iner¬ 
tial particle are dominated by viscous drag of the parent 
medium and is described by a force equation. For the 
sake of simplicity we disregard all other processes that 
may occur in the system (e.g. multi-phase transitions). 

The equation of motion for a particle takes the stan¬ 
dard form 

r = v, v = — (V — v), (1) 

T 

where r is the Stokes time of the particle, V(t) = 
V(r (t),t) is the velocity field induced by a blast at the 
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position of the particle. 

We assume that the flow field velocity V(r,f) can be 
approximately described by the Sedov-Taylor solution for 
a strong explosion [21] 


V=§*(C), (2) 

where £ = r/R(t) (and r = 0, t = 0 corresponds to the 
initial location and ignition time of the blast). 


R = 


(?f 


(3) 


where R(t) is the position of the shock front, E is the 
total energy released in the blast and p is the density of 
the medium. The dimensionless parameter /3 is function 
of the polytropic exponent 7 (for 7 = 7/5: /? ss 1.033 
[21] ) and the function T(£) can be closely approximated 
by its limiting value W(^) ~ I/7, ifO < £ < 1, and 
'k(C) = 0 otherwise [21] . The velocity of the shock front 
is given by the derivative of expression 




(4) 


For a particle located at an initial position r = ro away 
from the blast ignition point at t = 0, Eq. 0 reduces to 
the scalar form 


We can estimate r T from Eq. ([5]) in the short-time limit 
where the velocity term is insignificant: 


q r 

T (t + t 0 ) 


= 0. 


( 8 ) 


The solution of this equation can be represented in 
terms of Bessel functions [22], although the complete so¬ 
lution is cumbersome. In order to avoid dealing with 
expressions containing special functions we can deduce a 
simplified estimation of r T based on the following kine¬ 
matic consideration. 

After being hit by the front shock wave the particle be¬ 
gins to accelerate (driven by fluid drag), so its trajectory 
is given by the expression 


r = r 0 +-t , 


(9) 


where a is the particle acceleration. This expression re¬ 
flects that at r(t = 0) = ro and r(t = 0) = 0 (initially 
the particle is at rest). The acceleration term can be 
estimated by matching Eq. ([9]) with the analytical solu¬ 
tion of Eq. ^ or by directly substituting the expression 
0 into this equation (acceleration being the r term in 
Eq. @ evaluated at r = ro and t = 0), so a ~ q{ro/to)/ T - 
Then for t = r, Eq. § leads to a complete description 
of the reference position 


Tref = r 0 + r 0 


0. 
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( 10 ) 


which allows us to write the scaling law ([T]) in the follow¬ 
ing form 


r 

r H- 

T 


q r 

T (t + t 0 ) 


= 0 , 


(5) 


where q = 2/(57) and to is the time required for the 
shock wave to reach the particle, estimated from Eq. 0 



From Eq. ^ we can readily deduce a scaling law for 
particle displacement when inertial effects are negligible. 
By dropping the first term in this equation we arrive at 


r ft 

^ref \^ref 

where r re f is some arbitrary reference position which de¬ 
fines the particle location at time t = t re f. The r re f 
and t re f scales have been introduced to satisfy two ini¬ 
tial conditions of the original equation ([5]). If we define 
t re f = to + t ( i.e . time when inertial effects become 
unimportant) then r re f = ro + r T , where r T is the par¬ 
ticle displacement during the Stokes time r (i.e. from f 0 
to t 0 + t). 



r 

ro 


q t 


t 


to + r 


( 11 ) 


Eq. is the main result of the present study. We 
can see that the exponent in this scaling law of particle 
displacement depends only on the properties of the me¬ 
dia (since q = 2/(57)), an d is independent of both the 
properties of the particles and energy of the explosion. 
Moreover, since 7 > 1 |21j the particle dispersion is al¬ 
ways slower than the ballistic regime (i.e. q < 1). In 
general, the dispersion process can be characterized by 
two limiting cases, depending on the value of the ratio 
r/to in Eq. (fTT|) . For given characteristics of the explo¬ 
sion (energy E), particle and medium properties (Stokes 
time r), the ratio r/to can be associated with the initial 
position of the particle ro by introducing the scale 


r* 


= /3 


(wf 


( 12 ) 


For the particles initially located within the sphere r < 
r* (below we refer to this case as the ‘near field’), we 
arrive at a simplified form of the scaling law © 
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For the opposite case (r > r*, the ‘far held’) Q = 1 [23] . 
We remark that particle properties can influence the 
value of Q only in the near held region. 

At some point the particle motion described by 
Eqs. (13) will be terminated and the particle will come 
to rest. The time of this termination corresponds to a 
disappearance of the driving velocity V(r, t) in Eq. ([!]) , 
or a deviation of the function V(r, t) from the strong ex¬ 
plosion model ([2]) ( i.e . when the shock wave significantly 
dissipates). A simple estimation of this termination point 
can be deduced from the following arguments. 

It is well known that a spherical shock wave loses en¬ 
ergy and eventually transforms into a spherical acoustic 
wave GUES]. This transformation is governed by an in¬ 
terplay between the non-linear and dissipative processes. 
As linear acoustic waves cannot generate a persistent how 
[21] , it is apparent that a particle cannot be advected any 
further when this process begins to dominate. Assuming 
that the shock wave transformation is mostly due to non¬ 
linear effects, and applying the condition R = c to Eq. 
Q, we can readily deduce a stopping time 


t s = 


(If (If- 


(14) 


where c is the speed of sound in the media and all particle 
dynamics are confined to t <C t s . 

Analogously, one can deduce an estimation for t s when 
the shock wave transformation is given by a dissipation 
process MM- Scaling of the dissipation length is given 
by a diffusion law 8(t) ~ (i^f) 1 / 2 , where v* is a well- 
known aggregated dissipation coefficient determined by 
viscosity and thermal conductivity (see Ref. G3 §96). By 
equating R(t), Eq. to 5(t) we arrive at the following 
expression for the stopping time determined by dissipa¬ 
tion 


t 



(15) 


which yields a different t s value compared to Eq. (14). 
In the present study we assume the dissipation coeffi¬ 
cient v * is relatively small, and hence the stopping time 
is dominated by non-linear effects described by Eq. (14). 

Setting t=t s in Eqs. (13), we can derive the following 
scaling law for the maximum particle displacement 


^ a r p 0 E k r\ (16) 

ro 

with values p = —5/2, k = 1/2 + 2 /( 157 ), h = 1 — 2 /( 57 ) 
for the near field scaling and p = — 1 / 7 , k = 1 /( 37 ), ^ = 
0 for the scaling in the far field [25] . 

The scaling of particle displacement with Stokes time 
fmax CX T h provides insightful information on the effect 
of particle properties in the system, since for a spherical 
particle 


dpPp 

18 /i ' 


(17) 


TABLE I. Numerically recovered scaling exponents of physi¬ 
cal syst em parameters against relative particle displacement, 
Eq.(16 1 , for an explosion in air (7 = 7/5). 


Near Field (r < r*) Far Field (r > r*) 

Theory Recovered Theory Recovered 


p = -5/2 (-2.5) p = -2.47 p = -5/7 (-0.71) p = -0.72 
k = 25/42 (0.6) k = 0.58 k = 5/21 (0.24) k = 0.23 
h= 5/7 (0.71) h= 0.71 h = 0 (0) h = 0 


Here, d p is the diameter of the particle and p p p its 
density; p is the dynamic viscosity of the media. For 
instance, with all other parameters being equal, Eq.(16) 


predicts the particle displacement scales with the media 

' max (X p 


viscosity (in the near field) as r„— or // 1 + 2 /( 5 7). 


III. NUMERICAL RESULTS 


In order to validate our analytical predictions for the 
scaling laws of Eq. (16), we numerically solve Eq. ([ 5 ]) 
with parameter ranges of ro, E and r [26] . We use the 
stopping condition Eq. (14) to calculate a termination 


point and evaluate the relative particle displacements 
Tmax/f 0 and then estimate the scaling exponents p , k , h 
from the log-log plots. In each set of simulations we 
change only one parameter keeping all other parameters 
constant. The reported parameter values are selected 
to represent a large range of conditions which cross the 
near/far field boundary at r*, with enough data points 
to recover the predicted scaling exponents. Addition¬ 
ally, these values must reside within the time constraint 
to < t < i«. A multitude of parameter values recover the 
scaling laws, however we plot a single representational 
value to remove any possible ambiguity in the results. 

As a foundation we model explosions in a diatomic gas 
(air), for which 7 = 7/5 (21] and q = 2/7 in Eq. ([ 5 ]). The 
dynamic viscosity parameter from Eq. (17) is assumed 
to be a constant value of p ~ 1.983 x 10^® Pa/s (air at 
room temperature), and the particle density represents 
steel ball bearings p p = 7874 kg/m 3 for all reported re¬ 
sults. The results of analytical predictions and numerical 
simulations are summarized in Table [T] 

Figure [l] presents a typical output of our simulations. 
It depicts the scaling response of a particle’s initial po¬ 
sition ro as it is varied between 0.01 m and 5 m to the 
numerical solution of Eq. <§• Two regimes of particle 
dispersion (near and far field, see Eq. ©) are indicated 
via the solid and dashed lines respectively. Both results 
are in good agreement with the theoretical derivations 
presented in Table [I] 

Energy scaling can be determined in a similar manner, 
with all other parameters held constant, energy E is ad¬ 
justed from 452 kJ/kg to 45.2 TJ/kg, the result of which 
can be seen in Figure [2] Here again we recover two scales 
with good agreement to theory: the near field (solid line) 
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FIG. 1. Scaling dependence of the scaled maximum displace¬ 
ment rmax/ro of the numerical solution to Eq. |5]l against ro 
(o) with r* indicated as a vertical dotted line. Control pa¬ 
rameters for the presented data are E = 4.52 xlO 10 J/kg and 
t = 5.52xl0 -6 s. Recovery of a —5/2 dependence (solid line) 
in the near field and a —5/7 dependence (dashed line) in the 
far field are clearly visible. 


and far field (dashed line) against the numerical results. 

Finally, we investigate scaling with the Stokes time 
r in Figure [3j As r is a function of many parameters 
(see Eq. |17| ) , we fix all medium variables and alter only 
the particle diameter (also fixing particle density) over 
a 1 x 10~ 4 —10 pm range. Scaling in the near field is 
recovered with good agreement, and as predicted there 
is no dependency on r in the far field. 

To verify the scaling laws derived in Eq. (161 are gen¬ 
eral, we also compare results for noble gasses (ie. when 
7 = 5/3 and p ~ 1.956 x 10 -5 Pa/s for helium at room 
temperature). Here, the near field exponents are pre¬ 
dicted to be p = —5/2, k = 29/50, h = 19/25, which 
recover to p = —2.5, k = 0.57, h = 0.75 numerically. In 
the the far field, p = — 3/5, k = 1/5, h = 0, which also 
scale as expected recovering p = —0.62, k = 0.19, h = 0. 


FIG. 2. Scaling dependence of E against the scaled maximum 
displacement r maa; /ro of the numerical solution to Eq. @(o). 
Control parameters for the presented data are ro = 0.2 m and 
r = 5.52xl0~ 6 s. Recovery of a 25/42 dependence (solid line) 
in the near field and a 5/21 dependence (dashed line) in the 
far field are clearly visible. 



ln(r) 


IV. DISCUSSION 

The scaling laws established above can provide some 
predictions that can be validated experimentally. For in¬ 
stance, it has been observed that there exists a range 
of particle sizes for which the inertia of the particles 
combined with the decay of the blast wave allows them 
to overtake the primary shock front [T2[. This effect 
strongly depends on particle size and there exists a 
threshold (i.e. a particle size limit) below which this 
effect does not occur T2J. The proposed framework al¬ 
lows us to formulate a quantitative criteria for this phe- 


FIG. 3. Scaling dependence of r against the scaled maximum 
displacement r maa ,/ro of the numerical solution to Eq. © 
(o). Control parameters for the presented data are ro = 0.2 
m and E = 4.52 x 10 10 J/kg. Recovery of a 5/7 dependence 
(solid line) in the near field is clearly visible. Displacement is 
independent of r in the far field. 


nomenon. 

Consider the dispersion of particles by an explosion in 
air (7 = 7/5). The existence of an overtake event directly 
follows from a comparison of the scaling laws for the po¬ 
sition of the shock wave and particle displacement. The 
position of the shock wave scales as t 2 / 5 , Eq. ([ 3 ]), while 
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the particle displacement scales non-uniformly: initially 
(to -C r) it scales as t 2 with time, Eq. ([9| , and then slows 
down to oc t 2 / 7 at the large time limit, Eq. (13). This 


implies that at the large time limit the particle is always 
behind the shock wave and can only overtake it during 
the initial (inertial) stage. Since particle displacement 
during the initial stage is given by Eq. © , this leads to 
the following condition for the particle to be in front of 
the shock wave: 




/ t \ 2/5 

= ro(— J < r 0 + (a/2)(t-t 0 ) 2 , t > t 0 . 

(18) 

By introducing a new variable y = (t/to) 1 / 5 , this condi¬ 
tion can be recast to a non-dimensional form 


®(y) = y 2 - x(y 5 - i ) 2 - 1 < o, y> 1, 


(19) 


where x — (<?/2)(f q/t) <C 1 and $(1) = 0. 

The function $(?/) has a single real root y± satisfying 
the condition y > 1 with its approximate value y 1 ss 
(i/x ) 1/8 oc (t/ to) 1 / 8 (3> 1. This root determines the time 
when the particle catches up and ‘penetrates through’ the 
decelerating shock wave (since t = t\ =y\to). This time 
corresponds to the particle displacement rq =y\ro > ro, 
after which the particle decelerates and the shock wave 
overtakes it again. 

Similarly, the second time the shock wave overtakes 
the particle follows from Eqs. © and ([l3]) : 


n 


ft\ 2/5 ^ft-h\ 2/7 

U) 2ri+r " 0 (^r) ■ 


t > h, (20) 


or in a non-dimensional form (substituting the t\ =y\to 
and 7'i =y\ro parameters from the first crossing and y as 
defined above) 

Ay) = y 2 - Q(y 5 - y\f n - vl > o, v>vu (21) 

and T(yi) = 0. The positive real root of T(y) has an 
approximate value y 2 ~ Q 7//4 oc (r/f 0 ) 5/4 > 1. Both 
roots (2/1,2/2) are dependent on the characteristics of the 
explosion, and the parameters of medium and particle 
(via constants Q)- It is evident that the consistency 
condition, 7/2 > yi, always holds for sufficiently heavy 
particles and for particles initially located in proximity 
to the center of an explosion. 

In essence, we have identified three consecutive events 
where the particle and shock wave cross (one possible 
set of parameters which observes this phenomena is pre¬ 
sented in Figure |4]). At t = t 0 (r = r 0 ), the shock wave 
initially reaches a particle which is at rest. Driven by in¬ 
ertia, the particle overtakes the decelerating shock wave 
at t = ti > to (r = ri). Finally, &tt = t 2 > t± (r = r 2 ) 
the shock wave again catches up to the decelerating par¬ 
ticle and overtakes it. The zeroth-order estimates of t\ 


and t 2 derived from Eq. (19) and (21) respectively ob¬ 


tain an acceptable agreement with the numerical results 
of Figure [I] t x ss 2xl0 -8 (est) = 9xl0 -9 ( numeric) and 
t 2 « 2.5xl0~ 7 (est) = 9.5x10”' (numeric). 



FIG. 4. Shock wave trajectory (dashed) and particle position 
(solid) after an explosion with a yield E = 0.1 J/kg. A steel 
particle with diameter d p = 25 nm is initially at rest at posi¬ 
tion ro = 5 mm then starts accelerating at time to (o) after 
the shock passes over. The particle overtakes the shock at t\ 
(•), then at a later time t 2 (■) it is re-captured and remains 
behind the shock wave for the rest of the explosive event. 


An experimental validation of the ‘wave-particle’ over¬ 
take phenomena is a challenging undertaking which re¬ 
quires precise (and simultaneous) measurements of the 
positions of shock waves and tracer particles; as such 
there are few publications on this subject. To the best of 
our knowledge there is only a single experimental study, 
Ref- [121 in which the ‘wave-particle’ overtake phenom¬ 
ena has been observed and reported. The positions of 
the particles and shock in this study have been detected 
by means of two 150 kV flash X-ray pulsers and six piezo¬ 
electric pressure transducers respectively. Our interpre¬ 
tation of these phenomena, presented in this paper, stems 
from a simple kinematic analysis of the scaling laws for 
particle displacement, is in qualitative agreement with 
the experimental observations reported in Ref. [T2} 
Another interesting effect associated with an explosive 
energy release in a particle system is the formation of a 
residual sparseness (cavity) in an initially uniform parti¬ 
cle distribution ( e.g. dust) after the particles have been 
displaced by the shock wave from their initial position. 
We can easily estimate the scale of this cavity by invok¬ 
ing the analytical framework presented above. In fact, 
the edge of the cavity is formed from the far field parti¬ 
cles, for which to r. By setting t=t s and r ~ ro ~ r cav 
in Eq. (13) we can derive an estimate 


^ cav — 



( 22 ) 


where k = (4/3 5 /25) 1 / 3 ss 0.09, which can also be estab¬ 
lished based on dimensional arguments. This expression 
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provides a characteristic scale of the density of parti¬ 
cle distribution associated with the residual sparseness 
caused by an explosion (for a visual example of such a 
cavity, see images in Ref. ITUI and Figure 1 of Ref. ED). 


and defines the near field scaling law as 



V. CONCLUSIONS 

In summary, we have presented scaling laws for particle 
displacement caused by a rapid (and localized) energy re¬ 
lease (explosion). These scaling laws account for particle 
properties (mass, diameter, density), properties of the 
medium a particle is advected through (viscosity, den¬ 
sity, speed of sound) and the energy of the explosion dis¬ 
sipated through the shock front. We demonstrated that 
by employing a conceptual model of particle displacement 
that includes the strong explosion model and simplified 
particle kinematics, the exponents of these scaling laws 
can be derived analytically, which are in good agreement 
with numerical simulations. 

This framework has been constructed in a manner such 
that more realistic descriptions (different models of ex¬ 
plosion, non-Stokes drag, multi-phase transformations, 
etc ) can be implemented to ascertain improved estima¬ 
tions of each scaling law outlined here. For instance, tak¬ 
ing inertial effects of the fluid flow into account leads to 
an equation of motion with a different form of drag )21| 

r = v, v =i(V-v) 2 , (23) 

where the parameter A has a dimension of length and pro¬ 
vides an aggregated characteristic of particles and media 
(size, density, viscosity), similar to Eq. It is evident 
from this equation that the same scaling laws for particle 
displacement exist in the far field, viz ., with the expo¬ 
nents given in Table [IJ This and further additions will be 
elaborated upon in our future work. 

We anticipate that the results presented in this study 
will be useful in evaluating high fidelity models of particle 
transport and interpreting experimental observations of 
explosion phenomena. 


Appendix A: Near and Far Field Limits 


Eq. (11) is characterised by two limiting conditions, 
dependent on the ratio r/to- 

We initially investigate the ‘far field’, where t/ to <C 1 
and r > r* by taking the limit 


lim — 
r/to—tO r 0 



(Al) 


and defines the simple scaling law 

- = Q(rY> Q = L ( A2 ) 

r ° AW 

The ‘near field’ conditions are obtained via the condi¬ 
tion t/ t o 3> 1 (r < r»). This limit is given by 


lim — 

T/t 0 ->oo ro 



(A3) 


Appendix B: Scaling Law of Maximum Particle 
Displacement 


To derive the scaling law for maximum particle dis¬ 
placement in Eq. (16) we set t = t s in Eqs. (13) such 
that 



(Bl) 

(B2) 


= Q 




5/6 1 « 

- 5/3 - 6/2 
c ' o 


(B3) 


The general case, where q = 2/(5y), we arrive at 


rmax 

ro 


Q 




l 

37 


C 



(B4) 


For the far field, Q — 1 and we recover the scaling 
conditions 


rmax 

ro 


oc r n 7 E 3 tt u 


(B5) 


The near field requires further treatment as Q in this 
regime is non-scalar. 


Q = 


‘+1 


( 7 f) 5/2 (« I/2 . 


(?r<i) 


.1/2 


ro 

PJ 


n i)’ 

q fr 0 \ ~ 5/2+5q/2 ( p \ - 1 / 2 + 9/2 


v g/2 


— T 


2 

1/7 


(!) 


(B6) 


(B7) 


1/57 


(!) 


1 fr 0 \ 5/2+1/7 / p \ -1/2+1/57 _ x _ 

\e) 


57 V P 


-2/57 


(B8) 


Substituting Eq. (B8) into Eq. (B4), we find 
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T max 

ro 




'5 7 -I- 

57 



5+1. 

' -V 


2 

157 _ 2 

7" 57 


1 

57 





2 

37 


_ 5 

/3 3 ->' 



5 I 2_2_ _ 

/3 2_I “ 3 -TC 3 1T 0 


5 

2 


1- 

r 


2 

57 



1 _ 2 

2 157 


37 _ I 

c 3 ^r 0 7 


(B9) 

(BIO) 


This result reveals two scales in the near field, with a 
critical point at 


quoted in Table [T] of the main text: p = —5/2, k = 
25/42, h = 5/7. 


= 5q ^ro^ 5/2 ^P_y/ 2 = 57^0 


(Bll) 


representing the position at which the acceleration 
term dominates the initial position term of Eqs. 10 and 
11 in the main text. IV* is inversely proportional to r/to 
and as such IV* 1 represents very-small times, with 
IV* -C 1 being the primary term. 

We can now recover the scaling conditions 


Appendix C: Numerical Treatment of the Eq. © 
ODE 

The second order ODE of Eq. © was decoupled to a 
set of first order ODEs. Let yi = y and 2/2 = y, giving 
the first order system 


T max 

ro 


r^E 15 ~iT 


A7»l 


Tmax 

ro 


oc r n 


1 1 I 2 1 _2_ 

2 ]£ 2 "T" 157 q - 1 2 3 4 5 6 7 8 9 10 57 _ 


iV*<l 


(B12) 

(B13) 


m = 2/2 

1 

2/2 = --2/2 
T 


g 2 /i 

T (f+ t 0 )' 


(Cl) 

(C2) 


IV* 1 exists only in the limit t —> to and can therefore 
be dropped from consideration. 

Finally, as an example of scaling conditions, we can 
use 7 = 7/5 (air) to obtain the values for the near field 


This system was then solved in MATLAB using the 
ode 15s stiff solver, which implements a variable order 
(variable step size) method based on finite difference for¬ 
mulas. 
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